#pragma once

#include "bspline_opt_sep2_impl.h"
#include "bspline_opt_unsep_impl.h"
#include "bspline_opt_global_impl.h"
#include "controls_sample_impl.h"
#include "segmap_impl.h"

namespace outline {

/////////////////
// bspline_opt //
/////////////////

enum node_mobility {
	node_static = 0, ///< node with tangents is fully fixed
	node_cvar = 1, ///< center is variale
	node_lvar = 2, ///< left is variable
	node_rvar = 4, ///< right is variable
	node_indivisible = 8 ///< variables of the node can be optimized only together
						 ///< (either center is not fixed or node's type is smooth or symmetrical)
};

template <typename T, typename R>
void bspline_opt<T, R>::doit(
	segmap<T> const & target,
	std::vector<junc_type> const & junc,
	bspline_opt_params const & params,
	controls<R> & ctrls,
	std::vector<point> const * attractors,
	std::vector<R> const * tractions,
	std::vector<int> const * fixed_dir_tangents,
	std::vector<point> * poles,
	std::vector<point> * samples
)
{
	if (poles) 
		poles->clear();
	if (samples) 
		samples->clear();

	if (!attractors) {
		assert (!tractions);
		return core(
			target,	junc,
			params,
			ctrls, 0, 0, fixed_dir_tangents,
			poles, samples);
	}

	uint s = ctrls.nodes();
	bool closed = ctrls.closed();
	assert (junc.size() >= s);
	assert (closed && attractors->size() >= 3*s && tractions->size() >= 3*s ||
		    !closed && attractors->size() >= 3*s-2 && tractions->size() >= 3*s-2);
	assert (!fixed_dir_tangents || fixed_dir_tangents->size() == s);

	//determine type of each node: is it fixed or can it be optimized individually
	std::vector<int> mobility(s);
	bool all_static = true;
	bool all_indivisible = true;
	for (uint i = 0; i < s; ++i) {
		int & mob = mobility[i];
		mob = node_static;
		if ((*tractions)[3*i] < 1)
			mob |= node_cvar | node_indivisible;

		if (junc[i] != junc_flat) {
			if (i == 0) {
				if (closed && (*tractions)[3*s-1] < 1)
					mob |= node_lvar;
			}
			else if ((*tractions)[3*i-1] < 1)
				mob |= node_lvar;

			if ((closed || i != s-1) && (*tractions)[3*i+1] < 1)
				mob |= node_rvar;
		}

		if (is_smooth_or_symm(junc[i])) {
			int lr = (mob & (node_lvar | node_rvar)); 
			if (lr)
				mob |= node_indivisible;
		}
		all_static = all_static && mob == node_static;
		if (closed || i > 0 && i < s-1) 
			all_indivisible = all_indivisible && (mob & node_indivisible);
	}

	if (all_static)
		return; //nothing to do
	if (all_indivisible || s == 1)
		return core( //can't simplify task
			target,	junc,
			params,
			ctrls, attractors, tractions, fixed_dir_tangents,
			poles, samples);

	//the most complex case: nodes are subdivided on independent optimizable parts
	uint b, e;
	if (closed && ((mobility[0] & node_indivisible) || (mobility[s-1] & node_indivisible))) {
		b = s-1, e = 0;
		while ((mobility[b] & node_indivisible))
			--b;
		while ((mobility[e] & node_indivisible))
			++e;
		assert (b >= e);
		if (b == 0 && e == 0)
			return core(
				target,	junc,
				params,
				ctrls, attractors, tractions, fixed_dir_tangents,
				poles, samples);

		std::vector<junc_type> subjunc;
		subjunc.reserve(s - (b-e-1));
		subjunc.assign(junc.begin() + b, junc.begin() + s);
		subjunc.insert(subjunc.end(), junc.begin(), junc.begin() + e + 1);
		assert (subjunc.size() == s - (b-e-1));

		controls<R> subctrls;
		ctrls.get_subctrls(b, e, subctrls);

		std::vector<point> subattractors;
		subattractors.reserve(3*(s - b + e) + 1);
		subattractors.assign(attractors->begin() + b*3, attractors->begin() + s*3);
		subattractors.insert(subattractors.end(), attractors->begin(), attractors->begin() + e*3 + 1);
		assert (subattractors.size() == 3*(s - b + e) + 1);
		
		std::vector<R> subtractions;
		subtractions.reserve(3*(s - b + e) + 1);
		subtractions.assign(tractions->begin() + b*3, tractions->begin() + s*3);
		subtractions.insert(subtractions.end(), tractions->begin(), tractions->begin() + e*3 + 1);
		assert (subtractions.size() == 3*(s - b + e) + 1);

		std::vector<int> subfixed_dir_tangents;
		if (fixed_dir_tangents) {
			subfixed_dir_tangents.reserve(s - b + e + 1);
			subfixed_dir_tangents.assign(fixed_dir_tangents->begin() + b, fixed_dir_tangents->begin() + s);
			subfixed_dir_tangents.insert(subfixed_dir_tangents.end(), fixed_dir_tangents->begin(), fixed_dir_tangents->begin() + e + 1);
			assert (subfixed_dir_tangents.size() == s - b + e + 1);
		}

		core(target, subjunc,
			params,
			subctrls, &subattractors, &subtractions, 
			subfixed_dir_tangents.empty() ? 0 : &subfixed_dir_tangents,
			poles, samples);

		ctrls.set_subctrls(b, e, subctrls);
		b = e;
	}
	else
		b = 0;

	for (; b+1 < s; ) {
		assert (!(mobility[b] & node_indivisible));
		if (!(mobility[b] & node_rvar) && !(mobility[b+1] & (node_indivisible | node_lvar))) {
			++b;
			continue;
		}
		for (e = b+1; e < s && (mobility[e] & node_indivisible); ++e);
		if (e < s) {
			std::vector<junc_type> subjunc(junc.begin() + b, junc.begin() + e + 1);
			controls<R> subctrls;
			ctrls.get_subctrls(b, e, subctrls);
			std::vector<point> subattractors(attractors->begin() + 3*b, attractors->begin() + 3*e + 4);
			std::vector<R> subtractions(tractions->begin() + 3*b, tractions->begin() + 3*e + 4);
			std::vector<int> subfixed_dir_tangents;
			if (fixed_dir_tangents)
				subfixed_dir_tangents.assign(fixed_dir_tangents->begin() + b, fixed_dir_tangents->begin() + e + 1);

			core(target, subjunc,
				params,
				subctrls, &subattractors, &subtractions,
				subfixed_dir_tangents.empty() ? 0 : &subfixed_dir_tangents,
				poles, samples);
			ctrls.set_subctrls(b, e, subctrls);
		}
		b = e;
	}
}

template <typename T, typename R>
class bspline_opt<T, R>::append_samples_n_poles {
	segmap<T> const & target;
	std::vector<point> * poles;
	std::vector<point> * samples;
public:
	append_samples_n_poles(segmap<T> const & t,	std::vector<point> * p,	std::vector<point> * s) : 
		target(t), poles(p), samples(s) { }
	void operator () (uint, R, point const & pos, point const & tangent) const {
		point norm(-tangent.y, tangent.x);
		if (poles)
			poles->push_back(target.nearest(pos, norm));
		if (samples)
			samples->push_back(pos);
	}
};

template <typename T, typename R>
template <class S>
class bspline_opt<T, R>::append_eq {
	segmap<T> const & target;
	S & sp;
public:
	append_eq(segmap<T> const & t, S & s) : 
		target(t), sp(s) { }
	void operator () (uint piece, R t, point const & pos, point const & tangent) const {
		point norm(-tangent.y, tangent.x);
		sp.eq(piece, t, pos, target.nearest(pos, norm));
	}
};

template <typename T, typename R>
template <class C> 
void bspline_opt<T, R>::global_fit(
	segmap<T> const & target,
	std::vector<junc_type> const & junc,
	bspline_opt_params const & params,
	controls<R> & ctrls,
	std::vector<R> const * tractions
)
{
	bool closed = ctrls.closed();
	uint s = ctrls.nodes();

	//skip global fit if all nodes are fixed (but some tangents are not)
	if (tractions) {
		bool all_fixed = true;
		for (uint i = 0; i < s; ++i) 
			if ((*tractions)[3*i] < 1) {
				all_fixed = false;
				break;
			}
		if (all_fixed)
			return;
	}

	typedef gspline<C> gspline_t;
	gspline_t gsp(closed, s, junc, tractions);
	for (uint it = 0; it < params.iterations; ++it) {
		ctrls.sample(append_eq<gspline_t>(target, gsp), params.sampling_max_distance, !closed, true);
		gsp.eq_solve(ctrls);
	}
}

template <typename T, typename R>
void bspline_opt<T, R>::core(
	segmap<T> const & target,
	std::vector<junc_type> const & junc,
	bspline_opt_params const & params,
	controls<R> & ctrls,
	std::vector<point> const * attractors,
	std::vector<R> const * tractions,
	std::vector<int> const * fixed_dir_tangents,
	std::vector<point> * poles,
	std::vector<point> * samples
)
{
	bool closed = ctrls.closed();
	uint s = ctrls.nodes();

	///for the separable phase, give to the tangents with fixed direction traction equal to 1
	std::vector<R> stractions;
	std::vector<point> const * sattractors;
	if (tractions) {
		stractions = *tractions;
		sattractors = attractors;
	}
	else {
		stractions.resize(ctrls.size(), real(0));
		sattractors = &ctrls;
	}
	if (fixed_dir_tangents) {
		for (uint i = 0; i < s; ++i) {
			if (junc[i] == junc_flat)
				continue;
			bool left_dir_fixed = ((*fixed_dir_tangents)[i] & fixed_left_tan_dir) != 0;
			bool right_dir_fixed = ((*fixed_dir_tangents)[i] & fixed_right_tan_dir) != 0;
			if (is_smooth_or_symm(junc[i]) && (left_dir_fixed || right_dir_fixed))
				left_dir_fixed = right_dir_fixed = true;
			if (left_dir_fixed && i > 0)
				stractions[3*i-1] = 1;
			if (left_dir_fixed && closed && i == 0)
				stractions[3*s-1] = 1;
			if (right_dir_fixed && (closed || i+1 < s))
				stractions[3*i+1] = 1;
		}
	}

	bool require_unseparability = false;
	for (uint i = 0; i < s; ++i) {
		if (junc[i] == junc_flat)
			continue;
		bool left_dir_fixed = fixed_dir_tangents && ((*fixed_dir_tangents)[i] & fixed_left_tan_dir) != 0;
		bool right_dir_fixed = fixed_dir_tangents && ((*fixed_dir_tangents)[i] & fixed_right_tan_dir) != 0;
		if (junc[i] == junc_smooth)
			left_dir_fixed = right_dir_fixed = true; //fixed during unsep phase
		if (!tractions ||
			left_dir_fixed && i > 0 && (*tractions)[3*i-1] < 1 ||
			right_dir_fixed && (closed || i+1 < s) && (*tractions)[3*i+1] < 1 ||
			left_dir_fixed && closed && i == 0 && (*tractions)[3*s-1] < 1)
		{
			require_unseparability = true;
			break;
		}
	}

	if (params.fit == fit_global || params.fit == fit_global_then_local) {
		//Global motion
		switch (params.global_transform) {
		default:
			assert (false);
		case translation:
			global_fit< g2::basic_translation<real> >(target, junc, params, ctrls, tractions);
			break;
		case rigid:
			global_fit< g2::basic_rigid<real> >(target, junc, params, ctrls, tractions);
			break;
		case rigid_scale:
			global_fit< g2::basic_rigid_scale<real> >(target, junc, params, ctrls, tractions);
			break;
		case affine:
			global_fit< g2::basic_affine<real> >(target, junc, params, ctrls, tractions);
			break;
		}
	}

	if (params.fit == fit_local || params.fit == fit_global_then_local) {
		spline2 sp(closed, s, junc, ctrls, params.antiloop, sattractors, &stractions);

		if (sp.freevars() == 0) {
			assert (require_unseparability);
			uspline usp;
			usp.create(closed, s, junc, ctrls, attractors, tractions, fixed_dir_tangents);
			usp.update_smooth_nodes(ctrls, attractors, tractions);
			for (uint it = 0; it < params.iterations; ++it) {
				ctrls.sample(append_eq<uspline>(target, usp), params.sampling_max_distance, !closed, true);
				usp.eq_solve(ctrls);
			}
		}
		else if (require_unseparability) {
			uspline usp;
			usp.create(closed, s, junc, ctrls, attractors, tractions, fixed_dir_tangents);
			for (uint it = 0; it < params.iterations; it += 2) {
				if (it > 0)
					sp.update_constraints(junc, ctrls, params.antiloop, sattractors, &stractions);
				ctrls.sample(append_eq<spline2>(target, sp), params.sampling_max_distance, !closed, true);
				sp.eq_solve(ctrls, attractors, tractions);

				usp.update_smooth_nodes(ctrls, attractors, tractions);
				ctrls.sample(append_eq<uspline>(target, usp), params.sampling_max_distance, !closed, true);
				usp.eq_solve(ctrls);
			}
		}
		else {
			for (uint it = 0; it < params.iterations; ++it) {
				ctrls.sample(append_eq<spline2>(target, sp), params.sampling_max_distance, !closed, true);
				sp.eq_solve(ctrls, attractors, tractions);
			}
		}
	}

	if (poles || samples) {
		ctrls.sample(append_samples_n_poles(target, poles, samples), params.sampling_max_distance, !closed, true);
	}
}

} //namespace outline
